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Sharply  localized  distributions  of  whistler  waves  have 
been  observed  in  the  ionosphere  and  the  solar  wind.  In  this 
paper  it  is  demonstrated  that  such  distributions  could 
possibly  be  explained  by  a  modulational  instability  arising 
from  interactions  with  i on -eye l otr on  waves. 
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In  this  paper,  I  discuss  forced  integrable  systems  using 
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The  nonlinear  interactions  of  whistler  waves  with 
density  fluctuations  in  magnetized  plasmas  is  studied.  It  is 
found  that  this  will  be  a  major  interaction  in  proposed 
commercial  Tokamaks,  and  that  instead  of  non  1  i near i t i es 
dominating,  the  diffusion  effects  would  dominate  with  a 
subsequent  loss  in  the  coherence  of  any  such  whistler  beam. 
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Comment  on  “Specific  sine-Gordon  soliton  dynamics  in  the  presence 
of  external  driving  forces" 

D.  J.  Kaup 

Nonlinear  Studies  Institute.  Clarkson  College  of  Technology,  Potsdam,  New  York  13676 
(Received  27  September  1982) 

It  is  shown  that  the  results  presented  by  Reinisch  and  Fernandez  have  an  alternative  interpretation 
where  the  soliton  does  behave  as  a  Newtonian  particle.  The  key  features  required  for  this  alternative  inter¬ 
pretation  are  (i)  properly  defining  where  the  center  of  the  soliton  is,  and  (ii)  expanding  the  solution  so  as 
to  avoid  any  secular  terms.  When  these  two  objectives  are  achieved,  then  the  center  of  the  soliton  is 
found  to  satisfy  Newton's  equation  of  motion  for  a  point  particle. 


Recently  there  has  been  some  controversy  about  the  vali¬ 
dity  of  soliton  perturbation  theories  and  the  interpretation 
of  a  soliton  as  a  particle.  This  was  first  noted  in  Ko  and 
Kuehl’s1  study  of  the  K  dV  equation  with  time-dependent 
coefficients.  Their  result  for  the  position  of  a  soliton  xm, 
when  transformed  into  the  notation  of  Kaup  and  Newell,2 
gave 


dxm 

dt 


■w-f- 

3n 


where  2 19 2  is  the  amplitude  of  the  soliton  and  f  is  the 
damping.  On  the  other  hand,  a  soliton  perturbation  theory 
found2 

ft ■ 

where  jc  is  Kaup  and  Newell’s  position  for  the  soliton. 
These  sign  differences  are  real,  and  numerical  results  did 
support  Ko  and  Kuehl's  result.3 

More  recently  Reinisch  and  Fernandez  have  numerically 
studied4  the  sine-Gordon  kink  under  the  influence  of  a  con¬ 
stant  torque.  They  also  found  their  numerical  results  at 
variance  with  the  predictions  of  soliton  perturbation 
theories5*10  and  have  proposed  to  explain  this  by  declaring 
the  soliton  to  be  a  non-Newtonian  particle.  What  I  propose 
is  that  one  does  not  have  to  be  that  drastic  and  also  that  the 
Ko  and  Kuehl  observation  and  the  Reinisch  and  Fernandez 
observation  may  have  a  common  explanation. 

First,  let  me  state  some  facts,  then  I  shall  give  my  inter¬ 
pretation  of  these  results. 

(1)  The  soliton  or  kink  is  not  rigid  and  is  not  a  “point 
particle.”11  (Therefore  one  must  qualify  to  what  extent  one 
is  referring  to  it  as  a  “Newtonian”  or  a  “non-Newtonian” 
particle.  Should  one  look  at  the  short-time  or  long-time 
scales  to  see  this?) 

(2)  Any  “extended  panicle”  will  respond  with  a  time  de¬ 
lay  to  an  externally  applied  force.12  (This  is  also  verified  by 
Reinisch  and  Fernandez’s  numerical  results.  To  the  extent 
that  the  soliton  is  not  a  point  particle,  one  could  say  that  the 
soliton  was  non-Newtonian.  In  this  respect,  Reinisch  and 
Fernandez  were  correct.  What  they  observed  were  the  com¬ 
bined  transient  effects  of  a  soliton  reshaping  itself11  and  ex¬ 
periencing  a  time  delay.) 

(3)  The  expansion  used  by  Reinisch  and  Fernandez  con¬ 
tained  secular  terms.  (The  presence  of  secular  terms  limits 
the  validity  of  their  expansion  to  short-time  scales.) 

(4)  The  concept  of  a  “soliton”  comes  from  considering 


the  solution  for  t—  +«,  whereby  the  general  solution 
separates  into  “a  collection  of  solitons  in  a  sea  of  radia¬ 
tion.”  (Thus  to  identify  or  locate  a  soliton,  one  should  use 
a  solution  valid  for  large- time  scales,  not  short-time  scales.) 

(5)  The  definition  of  the  center  of  a  soliton  used  by 
Reinisch  and  Fernandez  is  different  from  the  definition 
used  in  soliton  perturbation  theories. 

Now  what  1  want  to  do  here  is  to  present  an  alternative 
interpretation  of  the  Reinisch  and  Fernandez  result.4  As 
they  did,  l  start  with  the  perturbed  sine-Gordon  equation 

un~  u„+sinu  —eR(x,t)  .  (1) 

I  now  difffer  from  their  procedure  and  instead  expand  u  as 

u(x.t)  “  £/0(x)  +  eulll(X,t)  +...  ,  (2) 

where  X  is  to  be  determined  and  (/0(x)  is  exactly  the  one- 
soliton  soliton  given  by 


U0(X)-4tan-1(e±l)  . 


(3) 


X  shall  be  defined  such  that  no  secular  terms  will  appear  in 
(2).  I  shall  define  the  center  of  the  soliton  to  be  at  the 
center  of  Uq ,  which  is  where  X  -0.  This  also  differs  from 
the  definition  of  the  center  according  to  Reinisch  and  Fer¬ 
nandez,  who  took  it  to  be  where  ti„(x,f)  was  a  maximum. 
To  avoid  relativistic  effects  and  to  maintain  simplicity,  I 
shall  take  X,  «0(e),  and  require  that 

X,2-l+X,2  .  (4) 

Then  from  (1),  (2),  and  (4),  the  first-order  result  is 

UotX,+fui"+*LuiU-(R  .  (5) 

where 

L  "  -d*+cost/o(x)  .  (6) 

The  operator  L  has  one  zero  eigenvalue,  which  is  a  bound 
state  whose  eigenfunction  is  proportional  to  £/0x-7  If  I  now 
demand  that  u<n  must  not  contain  any  secular  terms  in  this 
(first)  order,  then  «m)  must  be  orthogonal  to  this  bound- 
state  eigenfunction.  Thus  I  take 


u(l)(x,r)- J^dko»(t)/*(x)  , 
whence  both  X  and  ak  are  uniquely  determined  by 

r  Mx)R(x,t)dx 

|  **  " **  — 

£>)Uo,(x)dX 


(7) 


(8) 
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and 


which  evaluates  to 


J^mfk(X)R(x.t)dX  .  (9) 

In  the  above,  /»  and  /*  are  the  eigenfunctions  of  L  (Ref.  7) 
and  <u£-l  +k2.  Equations  (8)  and  (9)  are  the  results  for  a 
general  forcing  term  R(x,t).  As  one  may  see  from  (8),  the 
acceleration  X ,  of  the  soliton  is  directly  proportional  to  the 
bound-state  component  of  R(x.t),  while  from  (9),  the  am¬ 
plitude  of  the  continuous  spectrum  is  driven  by  the  fcth 
component  of  R(x,t). 

Now,  as  in  Ref.  4,  let  us  take  R  independent  of  x,  so  that 
we  have  a  constant  torque  being  applied  to  the  sine-Gordon 
field.  Also,  take  R  -0  if  t  <0  and  R  constant  for  t  >0. 
Then  for  t  >  0,  (8)  yields 

X.-  ±€-f*  ,  (10) 

4 


«ll,(x,/)-ir-*  2 - +  •  •• 


coshx 


(16) 


Now  from  (2),  (11),  and  (16),  we  have 


u(x.t)  -  Uo 


±-~tRt2\+^-R 


2- 


coshx 


+  0(t*)  . 
(17) 


Since  we  have  evaluated  (12)  by  a  Taylor’s  series  expansion 
in  z,  we  may  as  well  do  the  same  for  the  soliton  part,  noting 
that  l/o*“  ±2/coshX.  Whence 

u(x.t)-U0(x)+jtt2R+O(ti)  .  (18) 


from  which  we  obtain 

X«x  t^-tR-^-  .  (11) 

4  2 

Since  the  center  of  the  soliton  is  at  X -0,  it  then  follows 
that  the  soliton  (our  definition  of  the  center  at  least)  does 
behave  as  a  Newtonian  particle. JS'11 

However,  as  was  indeed  pointed  out  in  Ref.  4,  such  is  not 
observed.  And  to  understand  what  has  occurred  in  these  nu¬ 
merical  experiments,4  we  must  include  the  effects  of  the 
continuous  spectrum.  From  (9),  one  can  readily  obtain4'7 

u(l’(X.r)+-yfl  pjk  GU.x)[l-cosW>l  .  (12) 

where 


Jccos(*X)-sin(/lcx)ianhX 

G(k-X) — (r+P)«rr/i) 


(13) 


In  Ref.  4,  the  integral  in  (12)  is  evaluated  by  contour  in¬ 
tegration  and  is  reduced  to  an  infinite  series.  However,  that 
infinite  series  is  only  convergent  if  one  is  outside  the  light 
cone.  Inside  the  light  cone,  one  must  use  other  techniques. 
For  large  times  «m  will  conveniently  separate  into  the  two 
parts 

uU)--4rR  P  dkG(k,X) 

2  J-m 

+  fi?r  dk  G(k,  x)cos(mkt)  .  (14) 

2  */-« 


The  first  part  is  time  independent,  and  corresponds  to  a  per¬ 
manent  change  in  the  soliton's  shape.  The  second  term 
may  be  evaluated  by  stationary  phase,  and  represents  out¬ 
ward  traveling  radiation. 

If  instead  we  are  interested  in  short-time  scales,  then  we 
may  expand  (12)  in  a  Taylor  series,  obtaining 

-±Rt1fcmdkil+k1)G(k.x)  +0(t*)  ,  (15) 


Naturally,  Eq.  (18)  is  exactly  the  same  result  as  that  ob¬ 
tained  in  Ref.  4.  However,  1  have  obtained  it  via  a  different 
definition  and  interpretation.  I  interpret  Eq.  (17)  as  a 
Newtonian  particle  moving  with  a  constant  acceleration,  and 
with  radiation  being  created  at  a  rate  proportional  to  t2  on 
short-time  scales.  What  will  be  observed  numerically  is 
shown  in  Eq.  (18).  As  shown  by  Eq.  (18),  the  buildup  of 
the  above-created  radiation  will  exactly  cancel  the  soliton 
motion,  causing  the  soliton  to  seem  to  hang  motionless. 

One  can  also  explain  this  by  considering  the  various  time 
scale  involved.  In  a  Taylor  series  expansion  as  in  Eq.  (18), 
one  is  implicitly  considering  the  response  of  the  system  on  a 
very-short-time  scale,  at  least  faster  than  the  time  required 
for  a  signal  to  cross  the  width  of  the  soliton.  On  such  a 
short-time  scale,  for  example,  the  left  side  of  the  soliton 
will  not  know  what  has  happened  on  the  right  side  of  the 
soliton.  Thus  whatever  happens  on  the  right  side  cannot  ef¬ 
fect  the  left  side.  Thus  each  element  of  the  sine-Gordon 
field  will  respond  independently  of  all  other  elements.  To 
make  this  clearer,  consider  now  Eq.  ( 1 )  on  this  short-time 
scale,  and  in  the  rest  frame  of  the  soliton.  Since  we  start 
with  an  equilibrium  state,  we  have  ua  -  sin  u  -0,  at  least  on 
this  time  scale,  which  leaves  only 

u.-cRUr)  .  (19) 

What  (19)  demonstrates  is  simply  that  the  response  of  u  at 
x  is  independent  of  what  u  is  at  another  value  of  x  Each 
element  of  u  is  responding  like  a  free  particle,  independent 
of  ail  other  elements,  and  its  response  is  only  determined  by 
the  value  of  the  forcing  term  at  the  position  of  that  element. 
In  other  words,  the  concepts  of  solitons  and  radiation  are 
only  of  value  when  one  is  concerned  with  or  interested  in 
the  intermediate  or  long-time  behavior.  On  the  short-time 
scales,  the  soliton  concept  is  of  less  value  than  the  field 
concept,  as  was  demonstrated  by  Eq.  (19). 

I  also  suggest  that  a  similar  analysis  of  the  K  dV  equation 
may  well  explain  Ko  and  Kuehl's1  result,  but  that  remains 
to  be  seen. 
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It  is  found  that  nonlinear  coupling  of  whistler  waves  with  ion-cyclotron  fluctuations  can  lead  to  spiky 
whistler  electric  field  envelopes.  Analytical  results  for  the  latter,  which  are  solutions  with  discontinuous 
derivatives,  are  given.  Possible  application  of  our  results  in  solar-wind  plasmas  is  discussed. 


The  knowledge  of  nonlinear  propagation  of  electron 
whistler  waves  is  very  essential  to  the  understanding  of 
wave  phenomena  in  both  ionosphere  and  space  physics.1  In 
particular,  nonlinear  interaction  of  whistlers  with  adiabatic 
particle  motion,  ion-acoustic  oscillations,  and  magnetohy¬ 
drodynamic  perturbations  can  lead  to  the  modulational  insta¬ 
bilities. 3,3  The  latter  can  give  rise  to  self-focusing,  wave  lo¬ 
calization,  and  soliton  formation.3,4 

In  this  paper  we  consider  the  interaction  of  whistlers  with 
electrostatic  ion-cyclotron  oscillations.  It  is  found  that  such 
an  interaction  produces  spiky  whistler  envelope  solitions. 
Analytical  results  for  the  latter  are  given,  and  application  to 
space  plasma  is  discussed. 

Consider  the  propagation  of  right-handed  circularly  polar¬ 
ized  whistler  waves  in  the  form 

S-£(i  +/y  lexpOk  -  x-i<ui  )+c.c.  (1) 

The  frequency  <u  and  the  wave  number  k  are  related  by  the 
dispersion  relation 


where  no  is  the  unperturbed  density,  and  we  have  assumed 
oi/L  » <u(ft,-«).  Note  that  kill,  but  the  slow  modula¬ 
tion  is  in  both  x  and  z  directions.  Here,  Sn  and  v,  are  the 
electron  density  and  field  aligned  velocity  perturbations  as¬ 
sociated  with  the  low-frequency  plasma  motion.  The  two 
are  related  by  the  electron  continuity  equation 

9,8n  +  noSfU,  —  0  ,  (5) 

where,  in  view  of  the  low  frequency  (8, «  fl,),  we  have 
used  the  drift  approximation  vx -  vf -cr  x  ^  i<b/Bo, 
where  <b  is  the  slow  ambipolar  potential.  Note  that 
*9  ■  v£-0,  and  that  there  are  no  density  perturbations  as¬ 
sociated  with  the  whistler  waves. 

We  shall  limit  to  electrostatic  low-frequency  response,  and 
let  the  parallel  phase  velocity  of  the  modulation  be  much 
smaller  than  the  electron  thermal  velocity.  The  slow  elec¬ 
trons  are  then  in  equilibrium  along  B0.  From  the  parallel 
momentum  equation,  one  then  obtains  the  corresponding 
electron  density  perturbations3 


c3fc3  i  i _ _  (2)  8,8»/np—9,<h  +  F  ,  (6) 

oi  ^  6u  ( D  m  cosB  “  ) 

where  <P  —  e<b/T,,  T,  being  the  electron  temperature.  Here, 
where  cos •  -*,/*,  fc3-*3  +fc,3;  and  ft,  are  the  elec-  F  represents  the  low-frequency  ponderomotive  force  due  to 
tron  plasma  and  gyrofrequencies,  respectively.  The  cold  the  whistler  fields;  it  is  given  by3  5 
plasma  dispersion  relation  given  above  is  valid  for 
l<w  —  ft, I  »  k,vm,  in  »  ft i,  and  w  >  &>„/(l  +«aj»/fl,3)l/3,  ^ 

where  vm  is  the  electron  thermal  velocity;  ft,  and  <■»„  are  F  -  18, + 

the  ion  gyro  and  plasma  frequencies,  respectively.  '  * 


-La,  Igi2 
v,  16irn<>7", 


Nonlinear  interaction  of  whistler  waves  propagating  along 
the  external  magnetic  field  B0z  with  slow  background  plas¬ 
ma  motion  gives  rise  to  an  envelope  of  waves  governed  by 
the  nonlinear  Schrodinger  equation3,4 

i<d,  +  vfdI)£  +  4v;9,3£  +  77,8,3£-A£-0  ,  (3) 

where  vf  — 2 *>(ft,— «)/*ft,  is  the  group  velocity,  vj- v, 
x((lt-4<o)/k(l,,  and  T —  v,(ft,- 2«)/2iHO»  — «)  are 
the  parallel  (to  §o)  and  perpendicular  group  dispersions. 


Under  the  quasineutrality  condition,  the  slow  ion  motion 
is  given  by 

dfSn  +  ± ■  v(1  +  (iod,v.-0  ,  (8) 

8,v,j  -  -  c,1  V  ,d>  +  ft,va  x;  -  ^  i6n/n0  ,  (9) 

8,v»“ -e,3df$- vi9,8n/n«  .  (10) 

where  c, -  ( and  v, -  (y T,/m,)m  are  the  ion  sound 


respectively.  and  ion  thermal  velocities,  respectively.  The  effect  of  the 

The  nonlinear  frequency  shift  A  appearing  in  (3)  is  given  ponderomotive  force  on  the  ions  is  negligible  since  the 


by 


whistler  waves  do  not  involve  ion  motion. 


kv,  8 n  2v» 

2  n  o  v, 


Although  according  to  these  equations  the  ion  motion  is 
(4)  linear,  its  nonlinear  coupling  to  the  whistler  waves  occurs 
through  the  ambipolar  potential  <t>  appearing  in  (6).  From 
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(9)  we  get 

(9,2+ft2)va- -c,29,  i^-n.c,2  V  i<t> 

*z  —  ui9,  V  xSn/rto  .  (H) 

Combining  (8),  (10),  and  (11).  we  find 

192(92  +  ft,2)  —  vj V  29,2-  ft,2vi92]8n/no 

-c,2(  V29,2+n29I2)<t>  .  (12) 

where  V2-  V2+9/.  Eliminating  *  from  (6)  and  (12), 

we  obtain 

l92(92  +  ft  2)-c«  V  29  2-  ft  2c492l9,Srt/rt0 

--c,2(V29,2+ft,292)F  .  (13) 

where  c4  -  ( y  T,  +  T, )/ mi. 

Equation  (13)  describes  four  types  of  driven  low- 
frequency,  electrostatic  oscillations.  First,  for  9,-0  (adia¬ 
batic  response),4  we  obtain  from  (13)  that  v:-0,  and 

d,(&/i/n0)-lT./(T,  +  Tl)]F  .  (14) 

where  F  is  given  by  (7)  with  9,  -  0.  This  modulation  leads 
to  a  positive  density  perturbation  for  a ><  Cl,.  One  can  easi¬ 
ly  show  that  the  usual  hyperbolic-secant-profile  standing  en¬ 
velope  solitons  exist  if  a>  <  ft,/4  (i.e.,  vi  >  0). 

For  9,  «  ft,,  V  1  «  9,2,  and  ft,”2c42  V  J  «  1,  Eq. 
(13)  yields  the  driven  slow  ion-acoustic  oscillations,  name¬ 
ly,4 

(9,2-c429f2)9,S/i/»o- -cfflF  .  (15) 

For  this  case,  one  can  show  that  sub-  as  well  as  supersonic 
hyperbolic-secant-profile  solitons  exist.4 

For  ft  i  «  9,  and  9,2 «  V  one  obtains  the  driven  fast 
ion-acoustic  waves,  moving  predominantly  across  the  exter¬ 
na)  magnetic  field. 

In  this  paper,  we  are  interested  in  investigating  the  non¬ 
linear  state  involving  the  whistlers  and  the  ion-cyclotron 
fluctuations.  The  latter  satisfy  9, —  fl,  and  9r2 «  V 
Thus  we  obtain  from  (13) 

( 9 2  +  ft 2  -  ci V 2  )9,#n/no“  ~  c,29  ,  (16) 

which  gives  the  ion-cyclotron  wave  response  in  the  presence 
of  the  ponderomotive  force. 

We  look  for  stationary  solutions  of  the  coupled  set  (3), 
(5),  and  (16).  We  let 

E-W(f)expdSt)  ,  (17) 

where  W  and  8  are  real,  and  f  -  x  +  az  -  Vt.  The  constants 
a  and  y  are  real,  and  for  the  present  modulation  a  «  1. 

For  localized  perturbations,  we  require  v,(f ),  8 n({)>  and 
W(f)— 0  for  |f|— oo.  One  obtains  from  (5)  and  (16) 
u,”  v  fin  Inn  and 

8n/no“  ~Pt(  1  ~ 2vplv,)b\9l  .  (18) 

where  92  —  fV}/\6irnoT,,  v,-  Via ,  and  p,-c,/ft,.  Com¬ 
bining  (4)  and  (18),  we  get  A  -  Qdjf1.  where 

Q •»-j*v#(l -lVf/vgVpi  .  (19) 


The  wave  equation  (3)  then  yields  V  *au(,  and 

~Pd2t9-h9-Q9d\92-0  ,  (’0) 

where  P  —  a2vi+T.  Equation  (20),  which  contains  a 
derivative  in  the  nonlinear  term,  is  similar  to  that  describing 
upper-hybrid  solitons.11  It  has  a  localized  solution  given  by 
the  relation 

9  -  ffosechUx  +  az  -  ft  )/£  +  ( 1  -9l!9$  )1/21  .  (21) 

where  92  —  P/4Q  and  L2  —  P/ 28.  Note  that  8  >  0.  P  >  0, 
and  a  «  1  are  necessary  for  the  existence  of  the  soliton. 
Furthermore,  we  note  that  the  function  9  also  appears  in 
the  argument  of  the  hyperbolic-secant  function,  and  that  the 
theory  breaks  down  if  2v,  —  vt.  The  slope  of  the  soliton 
electric  field  9  can  be  shown  to  have  a  discontinuity  at  the 
center,  where  the  profile  consists  of  a  cusped  spike. 

Of  the  five  soliton  parameters,  namely,  the  amplitude  9 0. 
the  speed  V,  the  direction  cosine  a,  the  width  L,  and  the 
frequency  shift  8,  only  two  are  independent.  The  rest  are 
related  by  the  conditions 

2p}k  vg9 o  —  a2i/4  +  T  —  2L  28 

and  V  —  av,.  Here,  the  quantities  k,  vt,  and  T  are  charac¬ 
teristics  associated  with  the  whistler  carrier  waves. 

We  have  demonstrated  the  existence  of  sharply  peaked 
modulation  of  whistler  waves.  Such  modulations  originate 
from  the  nonlinear  coupling  of  whistler  waves  with  back¬ 
ground  plasma  motion  associated  with  electrostatic  ion- 
cyclotron  waves.  This  is  in  contrast  to  the  often-discussed 
whistler  modulation  by  ion-acoustic  or  magnetosonic  waves, 
where  the  resulting  pulses  have  smooth  hyperbolic-secant- 
profile  peaks.  The  physical  reason  for  the  discontinuous 
behavior  (of  the  derivative)  of  the  soliton  profile  may  be 
the  occurrence  of  a  marginal  balance  between  the  dispersive 
anil  nonlinear  effects,  in  the  sense  that  the  wave-breaking 
nonlinearity  is  barely  balanced  by  the  wave-spreading  disper¬ 
sion.  Similar  behavior  exists  in  certain  problems  in  hydro¬ 
dynamics.7  Physically,  the  spiked  profile  clearly  cannot  ex¬ 
ist.  It  is  expected  that  small  collisional  or  thermal  (Landau 
damping)  effects  will  smooth  out  the  profile  at  the  discon¬ 
tinuity  (of  the  derivative)  locally,  but  the  overall  peaked 
(spiky)  shape  will  remain. 

Discrete,  sharply  peaked,  whistler  wave  packets  seem  to 
have  been  observed  in  the  ionosphere*  and  in  the  upstream 
solar  wind.*”"  In  particular,  observations  in  the  foreshock 
region  of  the  Earth’s  bow  shock,"  as  well  as  near  inter¬ 
planetary  shocks  in  the  solar  wind,12  indicate  the  existence 
of  electron  or  ion  beam  related  whistler  turbulence,  which 
appears  as  spiky  pulses  of  electric  and  magnetic  field  fluc¬ 
tuations.  These  pulses  may  be  associated  with  the  ion- 
cyclotron  wave  modulated  whistler  wave  packets  investigat¬ 
ed  here.  The  lack  of  quantitative  data,  however,  prevents 
us  from  a  detailed  comparison  of  our  theory  with  the  obser¬ 
vations.  We  shall  nevertheless  compare  the  observed  pulse 
amplitude  with  our  calculations.  From  the  data  given  in 
Ref.  11,  namely,  no- 10  cm'\  7",**1  eV,  ck/a>  —  450. 

ft,/2»r- 140  Hz,  and  /-<u/2ir-56  Hz,  we  obtain 
Eilf  —  10'*  V2/m3Hz.  which  is  rather  high  in  comparison 
with  the  average  value  10”"  V2/mJHz  quoted  in  Ref.  II. 
The  difference  may  be  due  to  wave-dissipation  mechanisms 
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as  well  as  the  pulse  chain  behavior  (the  observed  pulses  are 
not  isolated  as  in  the  theory),  which  were  precluded  in  our 
theory.  The  condition  that  the  time  scale  of  the  modulation 
is  of  the  ion-cyclotron  frequency  ( •—  0.08  Hz  here)  is  con¬ 
sistent  with  the  estimated  (from  the  diagrams  of  Ref.  11) 
half-pulse-width  of  15  sec.  Unfortunately,  the  density  fluc¬ 
tuations  associated  with  the  pulses  do  not  seem  to  have 
been  recorded.  According  to  (18),  the  pulses  are  accom¬ 
panied  by  in-phase  density  humps  of  the  order 
8n//iQ—  10-4. 

Finally,  we  point  out  that  our  investigation  may  find  ap¬ 
plication  in  the  problem  of  radio-frequency  heating  of  plas¬ 


mas,  since  whistler  waves  (but  with  u>~  Cl,)  are  used  for 
electron  cyclotron  resonance  heating.  When  the  power  is 
sufficiently  high,  modulations  by  low-frequency  plasma 
motion  can  become  important.  The  ion-cyclotron  modula¬ 
tion  considered  here  would  be  of  particular  importance  if 
simultaneous  ion-cyclotron  resonance  heating  is  used  to 
heat  the  ions,  since  then  a  high  level  of  electrostatic  ion- 
cyclotron  fluctuations  would  be  present  to  enhance  the 
modulation. 
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A  method  for  solving  forced  integrable  systems  is  presented.  The  method  requires  the  knowledge 
of  at  least  one  piece  of  information  about  the  solution.  Once  this  is  known,  one  may  then  construct 
the  remainder  of  the  solution.  In  this  sense  these  systems  are  “almost  integrable.”  The  forced 
semi-infinite  Toda  lattice  is  used  as  an  example  and  to  illustrate  the  method. 

PACS  numbers:  03.40.Kf 


I.  INTRODUCTION 

Although  the  inverse  scattering  transform1  (1ST)  is  well 
established  as  a  method  for  solving  free  integrable  systems, 
little  work  has  been  done  on  forced  integrable  systems.  By 
“free”  we  mean  those  systems  without  some  type  of  forcing 
term.  Typical  examples  of  free  integrable  systems  would  be 
the  sine-Gordon  equation2 

<b„  -  <*x,  +  sin  <t>  =  0,  (1) 

where  the  boundary  conditions  are  <b  {x— ►  ±  ao  ,t )  =  2trn,  or 
the  nonlinear  Schrodinger  equation 

itl>,  ~ibxx±  2(i/>*ilr)if),  (2) 

with  the  boundary  conditions  of  tfix— ►  ±  ao,t)  — a,  where  a 
is  an  arbitrary  complex  constant.3,4  On  the  other  hand,  a 
“forced”  system  would  have  some  forcing  terms  which  de¬ 
termine  much  of  the  motion.  As  an  example  of  a  forced  inte¬ 
grable  system,  the  driven  sine-Gordon  chain  is  where  Eq.  ( 1 ) 
is  valid  for  x  >  0,  while  the  value  of  </>  (0,r )  is  externally  con¬ 
trolled.  If  one  drove  this  system  such  that  d>  (0,/ )  =  2m,  then 
for  every  one  unit  of  time  a  new  kink  would  have  been  inject¬ 
ed  into  this  sine-Gordon  chain.  Other  examples  are  easily 
imagined. 

One  will  note  that  the  above-mentioned  “free  integra¬ 
ble”  systems  are  all  completely  solved  by  the  1ST.  And  this 
method  of  solution  is  well  known.  But  in  general  the  “forced 
integrable”  systems  are  not  solvable,  except  in  special  cases 
wherein  one  may  utilize  some  symmetry.5  Otherwise  most  of 
what  we  know  of  such  forced  systems  has  been  obtained  by 
numerical  methods. 

If  one  reflects  on  what  happens  to  the  scattering  data  in 
a  forced  integrable  system,  one  can  appreciate  some  of  the 
complexity  of  such  systems.  For  example,  in  the  above-men¬ 
tioned  driven  sine-Gordon  chain,  the  scattering  data  must 
vary  as  some  complicated  function  of  time,  simply  because  in 
every  new  unit  of  time,  an  additional  kink  must  appear, 
which  means  that  a  new  pole  in  the  reflection  coefficient  has 
to  move  across  the  real  axis  up  into  the  upper  half  of  the 
complex  f-plane  (£  is  the  eigenvalue  of  the  scattering  prob¬ 
lem).  On  the  other  hand,  the  time  dependence  for  free  inte¬ 
grable  systems  is  quite  simple.  The  bound-state  eigenvalues 
are  fixed  in  time  as  is  also  the  magnitude  of  the  reflection 
coefficient.  Another  feature  of  these  forced  integrable  sys¬ 
tems  is  that  the  Lax  pair  relation6 


L,  =  [MX  ],  13) 

which  for  the  free  system  is  satisfied  everywhere  is  now  satis¬ 
fied  “almost  everywhere”  instead  of  “everywhere”.  Equa¬ 
tion  (3)  is  violated  at  those  points  where  the  system  is  being 
forced.  It  is  this  Lax  relation  which  guaranteed  the  integrabi- 
lity  of  the  free  system  in  the  first  place.  So  if  for  the  forced 
system  Eq.  (3)  is  now  satisfied  only  almost  everywhere,  could 
we  then  not  expect  such  systems  to  be  something  like  “al¬ 
most  integrable”?  Indeed,  such  is  the  case.  As  I  shall  demon¬ 
strate,  given  the  forcing  terms  and  only  a  few  additional 
pieces  of  information  about  the  system,  the  system  then  be¬ 
comes  completely  integrable.  This  additional  information  is 
not  independent  of  the  forcing  terms  and  is  quite  dependent 
on  them.  So  there  is  a  consistency  problem.  But  once  this 
additional  information  is  obtained  or  known,  then  the  re¬ 
mainder  of  the  system  does  become  completely  integrable. 

The  remainder  of  the  paper  will  be  devoted  to  using  the 
forced  Toda  lattice  as  an  example  of  an  almost  integrable 
system  and  to  illustrate  these  above  ideas.  By  “forced  Toda 
lattice”  I  mean  the  semi-infinite  Toda  lattice7 

i^a) 

(«>1), 

Pn=  -  exp (Q„  -  Q„  * , )  +  exp«2„  _ ,  -  CJ.  (4b| 

and  where  Q0  (and  P0  =  Q0)  are  externally  controlled.  In 
other  words,  Q0(t )  determines  how  the  zeroth  lattice  particle 
will  move  and  then  the  motion  of  all  other  particles  to  the 
right  of  this  particle  is  determined  by  Eq.  (4).  This  system 
was  suggested  to  me  by  Professor  Knopoff,8  who  along  with 
T.  G.  Hill9  had  observed  a  fascinatingly  regular  envelope 
structure  developing  out  of  an  apparently  chaotic  system. 
(See  their  Fig.  2.)  An  example  of  the  same  is  shown  in  my  Fig. 

1 ,  but  at  a  different  time.  What  one  should  note  is  the  regular 
envelope  structure  to  the  left,  whereas  as  one  moves  to  the 
right  the  structure  becomes  more  and  more  random  and 
chaotic.  To  say  the  least,  this  is  a  very  curious  and  strange 
behavior,  and  one  would  like  to  be  able  to  understand  what  is 
happening  here.  In  this  case,  the  forcing  of  the  zeroth  parti¬ 
cle  is  a  very  simple  uniform  forward  motion  Q0(t )  =  —  2b,  r, 
where  b,  is  some  negative  constant.  Thus  the  zeroth  particle 
is  being  rammed  into  the  other  panicles,  creating  a  shock 
wave.  The  strange  behavior  is  the  subsequent  creation  of  a 
regular  envelope  from  out  of  this  chaotic  shock  wave. 
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FIG.  1.  Plot  of*,  vjrin  the  forced  Toda  lattice  for  A,  =  —  1.95  at  f  =  64.0 
where  —  2*.  „  ,  is  the  velocity  of  the  nth  panicle.  Note  the  regular  envelope 
structure  to  the  right. 


The  study  of  shock  waves  in  one-dimensional  lattices  is 
not  new.  An  earlier  analysis  by  Holian  and  Straub 1 0  centered 
on  the  relaxation  toward  thermodynamic  equilibrium  in  the 
wake  of  shocks.  Included  in  their  numerical  analysis  was  the 
Toda  lattice.  These  numerical  results  for  the  Toda  lattice 
have  recently  been  intensively  analyzed1 1  by  using  local  1ST 
techniques.  (By  “local”  1ST  techniques  it  is  meant  that  one 
takes  a  small  section  of  the  system  and  analyzes  it  with  the 
1ST,  determining  what  solitons  are  present  inside  this  sec¬ 
tion,  etc.  Of  course  the  section  must  be  sufficiently  wide  so 
that  an  analysis  does  make  sense.)  This  is  in  contrast  to  what 
I  shall  do  here  which  would  best  be  described  as  a  “global” 
analysis.  Thus  my  analysis  is  a  compliment  to  theirs,  and 
many  of  our  results  are  of  course  the  same.  Mainly  we  differ 
in  emphasis.  Holian,  Flaschka,  and  McLaughlin11  sought  to 
explain  the  molecular-dynamics  experiments.  I  am  seeking  a 
more  general  method  for  determining  the  time  evolution  of 
the  scattering  data  when  an  integrable  system  is  being 
forced.  Only  the  model  and  the  specific  results  are  the  same. 
The  techniques  developed  by  each  of  us  are  different. 

Next  I  shall  briefly  summarize  the  1ST  for  the  semi¬ 
infinite  Toda  lattice  in  Sec.  II.  Then  in  Sec.  Ill  I  shall  deter¬ 
mine  the  time  dependence  of  the  scattering  data  for  the 
forced  Toda  lattice.  This  will  not  be  a  solution  of  the  initial- 
value  problem  since  this  solution  will  require  a  part  of  the 
solution  before  one  can  construct  the  problem.  So  there  will 
be  a  consistency  problem. 

Nevertheless  this  solution  is  still  useful,  and  in  Sec.  IV I 
shall  discuss  how  one  may  use  it  to  predict  the  scattering 
data  for  all  time.  I  shall  then  conclude  with  some  concluding 
remarks  on  the  consistency  problem. 

II.  THE  1ST  FOR  THE  FORCED  TODA  LATTICE 

Following  Flaschka, 12  we  define  a,  and  bH  by 

a,*.  =i«p[ -\{Qn (5a) 

I">11’  libl 


then  from  Eq.  (4)  it  follows  that 


an  =a,(h,  -bn_ ,). 
i>„  =2  (a2„  +  1  -  a2), 


(«>2), 


(6a) 

(6b) 


where  h,(/ )  and  Q0\i )  are  to  be  specified.  Equation  (6)  then 
determines  a,  and  b„  for  n>2. 

Consider  now  the  eigenvalue  problem12 

a,  +  tr„+1  +a„F_,  +  [bn-A  |F„  =0  («>1),  (7) 

where  A  is  the  eigenvalue  and  we  shall  take  a,  =  1  (see  Ref. 
13).  As  shown  by  Case13  one  may  define  the  scattering  data 
in  the  semi-infinite  discrete  case  as  follows.  (I  shall  shift  to 
the  AKNS  notation,  where  il>„  are  the  right  eigenstates  and 
are  the  left  eigenstates.)  Take 

A  =  J(z  +  l/z)  (8) 

and  assume  that  a„  —  1  and  bn  each  approach  zero  suffi¬ 
ciently  rapidly  that  the  following  results  hold.  Then  the  right 
eigenstate  may  be  defined  by 

i/>„(z)—z"  as  (9). 

where  ii>„z  ~ "  is  analytic  inside  the  unit  circle  of  the  2-plane.  I 
define 


Zjz^jl/z),  (10) 

which  is  the  second  independent  right  eigenstate  of  17). 

Now  define  a  left  eigenstate  by 

d„H=(2  -  l/zr'[  j/'o(2)l*„  (2)  -  I lr0(z)H>„  (Z)  ]  .  (11) 

By  construction, 

<*o  =  0,  (12a) 

*.-l.  (12b) 


Consider  using  Eqs.  (7)  and  (12)  to  construct  the  solution  <f>„ . 
Clearly  <*„  will  be  at  most  a  polynomial  in  A,  of  order  n  —  1. 
Thus  it  follows  that  d>„  is  analytic  in  A  except  for  a  finite- 
order  pole  at  A  =  oo . 

Define13 


S(z)  =  e2'«’'  =  0o(z)/ifio(z),  |13) 

where  5  is  the  phase  shift.  Then  the  scattering  data  consists 
of  the  values  of  6  (z)  for  z  on  the  unit  circle  (the  continuous 
spectrum)  and  the  poles  of  S  (z)  inside  the  unit  circle  (the 
bound-state  spectrum).  These  poles  are  the  zeros  of  \l>0{z)  in¬ 
side  the  unit  circle.  The  bound-state  part  of  the  spectrum  is 
specified  by  the  value  of  z  at  the  pole  (z,)  and  value  of  the 
normalization  constant  M],  which  is  the  negative  of  the  resi¬ 
due  of  z~'S  (z)  at  the  pole.  The  constant  M,  is  real,  whence 
M)>  0. 

The  inverse  scattering  equations  are  obtained  by  con¬ 
sidering  the  contour  integral 


7  2iri  ^o(z) 


l/zlz'"-1, 


(14) 


where  C  is  an  infinitestimal  circular  contour  CCW  around 
the  origin.  From  this  and  upon  expanding  tl>„  as 


<Uz)  =  AT„  £  Knjz>, 


(15) 


where  tm  =  l,  one  obtains  the  following. 1 3 1 2  First  construct 
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=T~<f-[t-SWI*.  <16* 

2m  J  z 

then  for  m  >  n>  1,  one  has 

*nm+Fn~m+  £  K^FJ  *  -  =  0  <17) 

j  —  n  +•  l 

from  which  one  may  solve  for  Next  construct  Km  from 

{K9r*  ml +FU  +  f  —  -  u«) 

Then  j,  and  fi„  may  be  recovered  from 

an=\K„/Kn.l,  1 19a) 

b,  =  J  (*„.„*,  (19b) 


From  these  equations  one  may  construct  the  direct  and 
inverse  scattering  transform  for  the  forced  Toda  lattice.  Giv¬ 
en  (a,  ,b„ )  for  n  >  2,  by  Eqs.  (7H 1 3 )  one  may  map  these  quan¬ 
tities  into  the  scattering  data.  And  given  the  scattering  data, 
from  Eqs.  ( 16H 19)  one  may  construct  the  inverse  scattering 
transform  which  allows  one  to  reconstruct  the  potentials 
(*„  ,bn )  for  2.  Clearly  we  may  do  either  of  these  at  any 
time.  Now  the  question  is,  if  ( a„ ,b„ )  for  n>2  evolves  accord¬ 
ing  to  Eq.  (6),  how  will  the  scattering  data  evolve?  This  we 
shall  answer  next. 


III.  THE  TIME  DEPENDENCE  OF  THE  SCATTERING 
DATA 

In  the  absence  of  forcing  and  when  one  has  an  infinite 
lattice,  Flaschlca'2  found  that  the  time  evolution  of  the  eigen¬ 
states  of  Eq.  (7)  was  given  by 

K -«.+  i  K„+1  -*„F„_,  +CV„,  (20) 

where  C  is  an  arbitrary  constant.  In  the  infinite  case,  the 
integrability  condition  for  (7)  and  (20)  is  the  infinite  Toda 
lattice  (Eq.  (6)  valid  for  all «].  But  in  the  semi-infinite  case, 
although  we  expect  Eq.  (20)  to  be  valid  for  large  n,  one  must 
carefully  account  for  the  equations  near  n  =  0  since  Eq.  (6)  is 
only  valid  for  n> 2.  Equation  (6)  just  cannot  be  true  for  n  =  1 
since  a,  and  bx  are  constrained.  Carefully  accounting  for 
these  equations  near  n  —  0  shows  that  for  the  forced  Toda 
lattice,  the  equivalent  form  of  (20)  is 

K  +  -*,K->+CK  (»> 2),  (21a) 

y^fC  +  A-bM-Vo,  (21b) 

K  =  (4*2  - 2i,)F,  +  v0{bx  -A  +  C).  (21c) 

We  comment  that  Eq.  (21b)  is  simply  Eq.  (21a)  for  n  —  1 
combined  with  Eq.  (7)  for  n  —  1.  Equation  (21c)  follows  upon 
differentiating  Eq.  (7)  with  respect  to  time.  One  may  easily 
verify  that  the  integrability  conditions  for  Eqs.  (7)  and  (21) 
are  now  Eqs.  (6). 

However,  one  may  not  uniquely  determine  the  time 
evolution  of  the  scattering  data  from  Eq.  (2 1 ).  Note  the  term 
present  in  Eq.  (21c).  From  Eq.  (5)  we  have 

a|  » i  exp ((?,  -  {?„),  (22a) 

fi,  *  -  2Go.  (22b) 

Although  we  do  know  h,  because  Q^t )  is  to  be  specified,  we 
do  not  know  what  * 2  will  be  because  £,(/ )  is  an  unknown. 


For  the  present,  let  us  assume  that  we  do  know  what  *2 
is,  and  continue.  To  determine  the  time  dependence  of  the 
scattering  data,  S  (z),  per  Eq.  ( 1 3)  we  require  the  time  depen¬ 
dence  of  ^2).  From  Eqs.  (9)  and  (2 1)  for  n  large,  I  determine 
that  for  the  eigenstate  (z),  the  constant  C  is 


C  =  —  ( |z  —  1/z).  (23) 

Define  the  function  \  \z,t )  by 

ti=Xe~  c,<  (24a) 

then  by  (21b), 

1>0  =  ec'[\A  -bx)x  -*]  (24b) 

and  (21a)  gives 

i+0-{z,  t)X~0.  (25) 

where 

a 1  -  4*2  -  i>,  -  |h,  -  A  )\  (26) 


Given  J2  :(z,/ )  and  the  initial  values  of  x  ir.O)  and  x  <*.0), 
one  may  construct  the  solution  for  x  ).  and  thereby  the 
solution  for  ri^z.t ).  From  Eq.  1 1 3)  one  may  now  construct  the 
scattering  function  5  (z,r ).  However,  the  value  of  4*2  (t  I  is  re¬ 
quired  before  any  of  this  may  be  performed.  If  4*2  U )  was 
known  then  the  remainder  of  the  solution  would  follow.  In 
this  sense,  these  forced  integrable  systems  are  “almost  inte¬ 
grate."  Some  piece  of  the  solution  must  be  provided  before 
the  remainder  of  the  solution  will  follow. 

However,  if  one  knows  something  nontrivial  about  the 
properties  of  4*2  (r ),  then  something  nontrivial  can  be  said 
about  the  scattering  data,  and  thereby  something  nontrivial 
about  the  remainder  of  the  solution.  It  is  in  this  manner  that 
I  shall  seek  to  glean  information  about  this  forced  system. 


IV.  THE  MOLECULAR-DYNAMICS  CASE 


Let  us  now  specialize  to  the  molecular  dynamics  case 
where  one  takes 


<?„(')  = 


0  if 
—  2b,t  if 


r<0, 

f>0, 


(27) 


with  fi,  as  a  constant,  —  2  bx  being  the  velocity  of  the  zeroth 
particle.  For  this  case  the  behavior  of  4*2  is  quite  simple1 1 
and  has  two  characteristic  forms.  These  are  shown  in  Fig.  2 
and  Fig.  3.  In  Fig.  2, 1  show  the  characteristic  form  of  4*2  for 
small  velocities;  in  this  example  fi,  =  —  The  main  features 
to  note  are  the  initial  rise,  followed  by  a  decaying  ringing, 
which  soon  decays  to  a  constant  value  of  approximately 
2.25.  The  value  of  bx  =  —  1  is  a  critical  value,"  and  for 
magnitudes  of  b ,  larger  than  this  critical  value  the  character¬ 
istic  form  of 4*2  changes,  as  one  can  see  in  Fig.  3.  Here  where 
h,  =s  —  2.0,  we  see  that  the  ringing  does  not  decay.  Instead 
4*2  seems  to  asymptotically  approach  an  oscillation  with  an 
amplitude  about  1.0  and  with  an  average  value  of  about  9.0. 

In  either  case,  the  dominant  feature  of  4*|  is  that  it 
shifts  from  1.0  at  /  =  0  up  to  some  larger  asymptotic  value, 
2.25  for  bt  =  —  0.5  and  9.0  for  b,  =  —  2.0.  So  as  a  first 
approximation  one  could  replace  4*2  in  Eq.  (26)  by  its 
asymptotic  average  value  and  then  proceed  to  solve  for  x 
from  Eq.  (25).  Of  course  this  will  not  generate  the  exact  solu¬ 
tion  for  the  scattering  data.  But  one  could  expect  that  it 
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FIG.  2.  A  plot  of  4d;  vs  t  when  0.5  showing  the  rapid  decay  of  the 

initial  ringing. 


would  contain  the  main  features  of  the  solution.  This  is  in¬ 
deed  so.  We  have  already  determined  that  this  procedure 
works  quite  well  for  predicting  the  soliton  birth  rate. 14 

As  a  final  point,  I  wish  to  point  out  that  there  may  be  a 
solution  to  the  consistency  problem  such  that  given  fi,(r ),  one 
may  be  able  to  directly  determine  4a; .  Let  me  illustrate  this 
in  the  molecular-dynamics  case,  Eq.  (27).  First,  I  determine 


the  initial  conditions  on  x  and  x-  At  t  =  0,  we  have 

a.  =J.  bn  =0  (n>2),  (28a) 

0i  -  J,  (28b) 

while  b,  is  some  nonzero  value.  Then  solving  (7)  for  gives 
tn  -z"  («>1).  (29a) 

to  =  1  -  2 b,z.  (29b) 

So  by  Eq.  (24)  we  have 


FIG.  3.  Aplotof*af  vs/ when  b,  -  -  2.0 showing  the  tsymptotic oscilla¬ 
tions. 


X(t  =  0)  =  z,  (30a) 

*(f  =  0)  =  }(r—  l)+zh,.  1 30b) 

Since  by  is  a  constant,  then  Eq.  (23)  may  easily  be  turned  into 
the  integral  equation 

2W  —  b,)x\t ) 

=  r  exp[(A  -  by\t )  +  (1  —  2b,z)exp[  —  (A  —  by)t ) 

+  2Jo  ^  ')sinh((A -h.Kr'-f)!,  (31) 

where  4a;  only  appears  in  the  kernel. 

Now  consider  the  analytical  properties  of  this  solution 
as  |r( — -0.  In  general  we  would  expect  essential  singularities 
at  z  =  0  due  to  the  presence  of  terms  like  e  £  But  now 
consider  (24a).  We  have 

(1/zltf,  =  [l/z\xea,  (32) 

where  ( 1  /z)ti  is  known  to  be  analytic  inside  the  unit  circle. 1 3 
For  arbitrary  values  of  4 a\  in  (3 1 ),  such  will  not  be  so  on  the 
right-hand  side.  One  may  easily  verify  this  by  using  a  Taylor 
series  expansion  about  t  =  0.  One  would  also  note  that  (32) 
would  have  the  correct  analytical  properties  only  if  4a3  satis¬ 
fies  the  equations  of  motion,  Eqs.  (6),  for  the  proper  value  of 
by.  (I  have  only  checked  this  out  to  second  order,  but  from  its 
form,  it  seems  reasonable  that  it  will  be  true  to  all  orders.) 

This  leads  us  to  conjecture  that  by  demanding  z~  'xec' 
to  be  analytic  inside  the  unit  circle,  the  correct  solution  for 
4a;  (f )  may  be  determined  and  obtained  without  having  to 
solve  the  equations  of  motion.  Given  fi,(r ),  Eqs.  (25)  and  (26) 
show  that  4oz  (r )  is  a  potential  for  x,  while  x  satisfies  a  Schro- 
dinger-like  equation  on  the  semi-infinite  interval  />0.  Clear¬ 
ly,  Aa\  could  be  mapped  into  the  scattering  data  for  the  prob¬ 
lem  given  by  Eq.  (25).  But  whether  or  not  the  required 
analytical  properties  of  ^  in  Eq.  (32)  are  sufficient  to  obtain 
this  scattering  data  remains  to  be  seen. 
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The  soliton  birth  rate  in  the  semi-infinite  Toda  lattice  is  studied.  The  lattice  is  forced  by  driving 
the  zeroth  panicle  with  a  constant  velocity  into  the  remainder  of  the  lattice.  An  approximate 
solution  for  the  soliton  binh  rate  is  derived  and  it  is  shown  to  compare  quite  favorably  with  the 
actual  binh  rate. 

PACS  numbers:  03.40. Kf 


In  a  recent  paper1  one  of  the  authors  (DJK)  discussed 
and  demonstrated  how  one  could  solve  “almost  integrable” 
systems,  one  example  of  which  is  the  forced  Toda  lattice. 
This  is  the  semi-infinite  Toda  lattice2  where  the  equation  of 
motion  is 

Qn  =  «P<G,  -  Qn - 1 )  -  exp(& * ,  -  {?„)  (1) 

for  n>  1.  The  position  of  the  zeroth  panicle  Q0{t )  is  assumed 
to  be  driven  by  some  external  agent.  And  the  motion  of  this 
particle  then  drives  all  other  panicles  through  Eq.  ( 1 ).  A 
simple  example  is  a  case  from  molecular  dynamics3  where 
one  starts  with  a  static  lattice;  then  at  t  =  0  one  forces  the 
zeroth  panicle  to  ram  into  the  remainder  of  the  lattice  by 
imposing  upon  it  a  uniform  forward  velocity.  Thus  Q0  =  v<yt, 
where  v0  is  a  constant. 

As  this  zeroth  panicle  rams  into  the  remainder  of  the 
lattice,  a  shock  wave  is  created,  the  front  part  of  which  con¬ 
sists  of  a  collection  of  solitons,  all  with  approximately  the 
same  velocity.  Parts  of  this  shock  wave  have  been  analyzed3 
using  “local  1ST”  techniques  to  verify  that  solitons  are  pres¬ 
ent  with  approximately  the  same  velocities. 

With  the  recently  developed  method  for  handling  al¬ 
most  integrable  systems, 1  it  now  becomes  possible  to  accura¬ 
tely  predict  what  the  soliton  structure  and  spectrum  of  this 
shock  wave  is.  The  purpose  of  this  paper  is  to  predict  this 
spectrum  and  to  compare  the  predicted  soliton  spectrum 
with  the  actual  observed  spectrum.  As  we  shall  see,  the 
agreement  between  the  predictions  and  the  numerical  results 
is  quite  good  indeed. 

Next  we  shall  summarize  those  equations  and  results 
from  Ref.  1  which  are  applicable  to  the  motion  of  the  soliton 
spectrum.  The  solution  of  these  equations  requires  one  to 
know  beforehand  what  will  be  the  separation  between  the 
first  two  panicles  as  a  function  of  time.  We  approximate  this 
in  a  reasonable  manner  and  obtain  thereby  an  approximate 
solution  for  the  motion  of  the  bound-state  (soliton)  spec¬ 
trum.  We  next  numerically  compute  the  lattice  motion  from 
Eq.  (1),  determine  what  the  actual  spectrum  is  at  various 
times,  and  then  compare  results. 

According  to  Kaup,'  the  inverse  scattering  transform 
(1ST)  for  the  forced  Toda  lattice  requires  the  solution  of  the 
eigenvalue  problem 

a,*,  6,+  ,  +  «.  #,-■  +(*„ ->1)^  =0,  (2) 

where 


A  =  \{2  +  1/Z),  (3) 

and  d>n  is  the  eigensolution  where 

\t>n  — *  z"  as  n  — »  -f  ao .  (4) 

The  quantities  an  and  b„  in  Eq.  (2)  are  related  to  Q„  by 

=}«p[  -  }((?„ -C?„-i)]>  (5a) 

(5b) 

and  thus  as  n  — ►  +  oo , 

—  J,  (6a) 

0.  (6b) 


Note  that  Z>,  is  just  the  negative  of  one  half  of  the  velocity  of 
the  driven  zeroth  panicle.  Also  a,  cannot  be  defined  by  (Sa) 
since  the  n  =  —  1  panicle  does  not  exist.  Instead  we  may 
define  it  to  be  J,  as  was  shown  by  Case.4  The  bound-state 
eigenvalues  are  those  values  of  z  where  #yz|  is  zero.4-5  These 


only  occur  when  z  is  real  and  is  between  —  1  and  +  1. 

As  shown  by  Kaup,1  if  one  defines  the  function  x  hy 

(7) 

it  then  follows  that 

X-U-b,)x-^~a,  (8) 

and  that  x  will  satisfy 

x+n2x= o,  (9) 

where 

/2 2(z,r )  =  4o2(r )- (h, -A)2- 6,  (10) 

and 

C=  —  |(z  —  l/z).  (11) 


Thus  if  one  possessed  the  function  x  (z,r ),  from  (8)  one  could 
construct  t )  thereby  obtaining  the  soliton  spectrum  (the 
zeros  of  tf>0)  as  a  function  of  time.  However,  before  we  may 
construct  the  solution  for  x>  we  must  know  4a\,  which  by 
(Sa)is 

4flJ=exp (Co-G.l-  (12) 

Although  Q0  is  given,  Q ,  is  not  and  requires  the  solution  of 
the  problem  which  we  are  trying  to  solve.  For  the  moment 
we  shall  simply  assume  that  4 a\  is  known,  and  continue. 

Assuming  4a\  (t )  to  be  known,  then  we  may  solve  Eq.  (9) 
as  follows.  Take  a  solution  of  (9)  to  be  of  the  form 
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X  =  A***\  (13) 

then  Eq.  (9)  gives 

A  =  const/)  <p  ),/2,  (14) 

d2  =  /}2+M2+/i.  (15) 

where 

M=A/A.  (16) 


The  initial  values  for  \  and  X  follow  from  the  initial 
values  for  a„  and  b„  as  follows.  Consider  the  moment  just 
after  t  —  0  where  h,  has  reached  its  nonzero  uniform  value. 
Here 


a.  =i  («>1),  (17a) 

b,  *0  (n>2).  (17b) 

We  may  now  solve  (2)  for  the  initial  value  of  \l>„ .  We  find 


4(r  =  0)  =  z*  («>1),  (18a) 

4(*  =  0)  =  1  —  2h,  z,  (18b) 

which  by  (7)  and  (8)  give  the  initial  values 

*  (f  =*=  0)  =  z,  (19a) 

X{t*=0)*A2-  1  +  h,z.  (19b) 

Matching  the  two  possible  solutions  in  (13)  to  these  initial 
conditions,  we  determine  the  correct  solution  of  x  to  be 

,  =  -Mzcosd  +  Xo-E.^.  sin  6 1,  (20) 

Aq  *■  00  J 


where  the  subscripts  “0”  refer  to  initial  values  and  we  have 
taken 


4*0.  (21) 

So  far  no  approximations  have  been  made.  From  (8),  the 
zeros  of  4  will  be  where 

*M-b,  ~fi)  +  i(l/z-A-b,  +  /*„) 

(/l  —  h,  —  ^)(l/z-yl  -by  +  fi0)  -4>o4>  ; 

Define 

d0  =  arc  tan  [  4/(  1/z  -A  -b,  +  fi0)}  (23) 

with  which  Eq.  (22)  can  be  reduced  to 

d  =  d0  +  arctan{d44  —  bx  —  n)].  (24) 

Now,  let  us  approximate  in  the  spirit  of  the  WKB  meth¬ 
od  to  determine  d  and  4-  From  ( 1),  (3),  ( 10),  ( 14),  and  (16)  one 
has  that  the  initial  value  of  n  is 

/i0-6,/2dg.  (25) 

Provided  4  was  not  close  to  zero,  the  solution  of  ( 15)  would 
be  ip  =  ±fi.  However,  if/}2  would  be  close  to  zero  we 
would  have  to  account  for  the  terms  fi1  +  ft.  We  do  this  by 
evaluating  them  for  4  small.  Otherwise  they  would  have  no 
significant  effect  and  could  be  ignored.  For  small  4,  we  have 

^=*hJ/*S«4^  (26) 

so  we  approximate  Eq.  (IS)  initially  by 

+  (27) 

which  is  a  cubic  equation  fordo-  It  has  only  one  positive  real 
root  when  A  and  b,  are  real. 

For  the  later  times,  we  shall  simply  ignore  the  effects  of 
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FIG.  I.  A  plot  of  4a:  vs  /  when  b,  =  —  0.5  showing  the  rapid  decay  of  the 
initial  ringing. 


H  and  fi.  Thus  we  take  them  to  be  zero  in  ( 1 5)  and  (24).  It  only 
remains  to  specify  the  values  of  4a;  (f ).  To  see  how  best  to  do 
this  consider  4a^  vs  t  as  shown  in  Fig.  1,  which  is  when 
b{  =  —  0.5,  and  Fig.  2  which  is  when  b,  =  —  2.0.  What  we 
observe  there  is  that  4a^  rapidly  shifts  from  its  value  of  1 .0  at 
t  =  0  to  a  larger  average  value.  Clearly  the  most  dominant 
feature  is  this  definite  shift  in  the  average  value.  So  we  shall 
approximate  the  value  of  4a2  required  in  the  calculation  of  d, 
Eq.  ( 1 5),  by  its  average  value.  Thus 

d2^<4a2>-(h,-/i)2.  (28) 

From  Figs.  1  and  2,  we  have 

6,  =  —  0.5,  <4a2>~2.25,  (29a) 

b,  =  -2.0,  <4a2>~9.0,  (29b) 
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FIG.  3.  The  solium  birth  rate  when  6,  —  —0.3  as  predicted  by  Eq.  (311 
(solid  linel  and  as  actually  is  (dashed  line). 


which  are  the  only  values  that  we  shall  consider  here. 

Now 

6  =  it,  (30) 

and  Eq.  (24)  gives 

f  =  ( ^ )~  *[.d0  +  arctan(  i  /[A.  —  h,))J.  (31) 

This  equation  gives  all  possible  times  associated  with  a  given 
possible  value  for  a  bound-state  eigenvalue. 

A  plot  of  these  t  values  vs  z  is  shown  by  the  solid  lines  in 
Fig.  3  for  h,  =  —  0.5  and  in  Fig.  4  for  h,  =  —  2.0.  In  Fig.  3, 
these  curves  are  easily  interpreted  as  being  the  motion  of  the 
eigenvalues  of  individual  solitons.  The  first  soliton  is  created 
at  r— 0.1  with  an  eigenvalue  ofjustabovez  =  —  l.fAsoliton 
with  z  —  —  1  would  have  a  zero  velocity,  zero  amplitude, 
and  an  infinite  width.  When  z  is  just  greater  than  —  1,  then 
these  values  become  finite  and  nonzero.)  This  eigenvalue 
moves  rapidly  toward  the  limiting  value  of  —  0.29  at  which 
all  bound-state  eigenvalues  eventually  tend  to  collect,  as  seen 
in  Fig.  3.  The  motion  in  Fig.  4  is  quite  similar,  except  that  the 
solitons  are  created  at  a  faster  rate,  and  the  first  soliton  al¬ 
ready  exists  at  t  =  0.  The  limiting  value  is  now  —  0. 10, 
which  means  faster  and  narrower  solitons  as  one  would  ex¬ 
pect. 

To  see  how  good  these  predictions  are,  let  us  compare 
this  with  the  actual  soliton  spectrum.  To  determine  this,  we 
shall  numerically  integrate  Eq.  (1)  up  to  some  time  t.  At  this 
time,  we  shall  calculate  the  a,  ’s  and  the  bH ’s  as  given  by  Eq. 
(5),  then  solve  Eq.  (2)  numerically  for  tf0(z),  plotting  i>^z)  vs  z 
from  z  =  —  1  to  z  —  +1.  One  may  then  easily  pick  out  the 
zeros  of  tfijz)  which  are  the  bound-state  eigenvalues. 

The  result  of  this  are  the  dashed  lines  in  Figs.  3  and  4. 
As  seen  in  Fig.  3,  the  agreement  is  quite  good,  the  only  differ¬ 
ence  being  slight  phase  shift  in  the  initial  birth  times.  Other¬ 
wise  the  eigenvalue  motion  is  quite  accurately  predicted  by 
Eq.  (24).  Figure  4  does  not  show  as  good  an  agreement,  al¬ 


FIG.  4.  The  soliton  birth  rate  when  b,  m  —  2.0  as  predicted  by  Eq.  (311 
(solid  linel  and  as  actually  ts  (dashed  linel. 


though  the  general  shape  and  motion  still  quite  accurately 
reflect  the  actual  curves.  This  discrepancy  may  arise  in  part 
from  ignoring  the  oscillations  in  4a\  (see  Fig.  2)  which  do  not 
seem  to  be  decaying  away.  They  do  rapidly  decay  away  in 
Fig.  I,  and  for  that  value  of  the  results  shown  in  Fig.  3 
gave  excellent  results. 

In  conclusion  we  have  demonstrated  that  one  can  solve 
for  the  soliton  spectrum,  and  its  subsequent  motion,  when  an 
integrable  system  is  driven  by  forcing  terms.  The  method 
does  require  having  some  particular  information  about  the 
solution,  so  it  is  not  a  method  for  solving  the  initial-value 
problem.  However,  the  information  required  for  finding  the 
soliton  spectrum  need  not  be  detailed,  and  we  found  average 
values  to  be  adequate  to  reproduce  at  least  the  gross  features 
of  the  curves. 
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Whistler  scattering  from  density  fluctuations  in  magnetized  plasmas 
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Scattering  of  a  coherent  whistler  from  random  density  fluctuations  is  treated  by  a  perturbation 
procedure.  The  attenuation  length  is  calculated  for  scattering  by  a  drift-wave  type  fluctuation 
model  and  is  simply  expressed  as  a  function  of  the  propagation  angle  in  the  limit  where  the 
whistler  wavelength  is  long  compared  to  a  correlation  length.  For  typical  tokamak  and  space 
plasmas,  this  scattering  becomes  most  important  at  larger  angles. 


I.  INTRODUCTION 

Wave-scattering  from  random  fluctuations  in  a  plasma 
has  been  a  subject  of  extensive  investigation  since  fluctu¬ 
ations  of  various  sorts  such  as  in  the  density,  temperature, 
magnetic  field,  etc.,  are  always  present  in  a  plasma.  For  in¬ 
stance,  electrons  in  the  upper  atmosphere  are  distributed  in 
an  irregular  fashion.1  Recently  density  fluctuations  have 
also  been  detected  in  a  tokamak  type  plasma.2'3  Thus  it  be¬ 
comes  essential  to  investigate  and  understand  the  influence 
of  random  fluctuations  on  wave  propagation. 

An  early  attempt  to  consider  the  effect  of  ionospheric 
irregularities  on  wave  propagation  was  due  to  Budden.4  He 
assumed  a  Gaussian  model  for  the  irregularities  and  suggest¬ 
ed  that  the  relatively  few  ground  observations  of  whistlers 
could  be  explained  in  terms  of  the  scattering  from  small  irre¬ 
gularities  having  a  finite  elongation  along  the  geomagnetic 
field.  Radio  wave  scattering  from  ionospheric  irregularities 
was  considered  by  Simonich  and  Ye.3  They  considered  a 
"bistatic”  problem  in  which  a  radio  wave  is  sent  out  into  the 
upper  atmosphere  by  a  transmitter  and  the  scattered  signal  is 
received  back  on  earth.  They  could  then  regard  the  region  of 
scattering  as  small  compared  to  the  distance  traveled  by  the 
radio  wave.  Thus  they  used  the  single  scattering  or  the  first 
Born  approximation  to  calculate  the  scattered  power  and  the 
scattering  cross  section  as  a  function  of  the  random  inhomo¬ 
geneity  of  the  medium.  Ott6  considered  the  scattering  of  a 
lower-hybrid  wave  from  random  density  fluctuations  using  a 
wave  kinetic  formalism.  Using  the  same  formalism,  Hui  et 
at.1  reported  a  numerical  calculation  of  the  ray  trajectory 
modification  of  the  electron-cyclotron  wave  due  to  the  den¬ 
sity  fluctuations. 

In  this  work  we  shall  consider  the  scattering  of  a  coher¬ 
ent  whistler  wave  from  random  density  fluctuations  in  a 
magnetized  plasma.  We  are  principally  motivated  by  the  fol¬ 
lowing  considerations:  (i)  Whistlers  are  known  to  travel  long 
distances  in  the  upper  atmosphere  being  guided  by  the  field 
aligned  columns  of  ionization.*  Since  the  electron  distribu¬ 
tion  in  the  upper  atmosphere  is  known  to  be  irregular, 1  whis¬ 
tler  propagation  in  the  presence  of  these  irregularities  would 
be  an  important  study,  (ii)  In  the  supplementary  radio  fre¬ 
quency  heating  scheme  in  the  lower-hybrid  range  of  frequen¬ 
cies.  it  has  been  suggested  that  the  use  of  the  "fast”  or  whis¬ 
tler  wave  may  have  some  advantages  over  the  usually 
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considered  “slow”  or  lower-hybrid  wave.910  Since  the  den¬ 
sity  fluctuations  also  occur  in  tokamaks,2-3  it  is  then  essential 
to  consider  whistler  propagation  in  the  presence  of  fluctu¬ 
ations. 

We  shall  follow  a  perturbation  procedure  due  to  Kel¬ 
ler'  1  to  study  the  whistler  scattering  problem.  An  important 
advantage  of  this  method  is  that  one  can  evaluate  the  final 
results  for  both  the  long  as  well  as  the  short  wavelength  lim¬ 
its.  Treating  the  fluctuations  to  be  weak  on  a  homogeneous 
background  plasma,  first  a  modified  dispersion  relation  in 
the  presence  of  fluctuations  will  be  derived.  It  will  be  seen 
that  the  fluctuations  cause  the  wave  refractive  index  to  be 
complex,  with  the  real  pan  representing  the  change  in  the 
phase  velocity  of  the  wave  and  an  imaginary  pan  signifying 
attenuation  of  the  primary  wave  due  to  scattering  from  fluc¬ 
tuations. 

Keller’s  method  was  employed  by  Liu'2  to  study  the 
scattering  of  the  ordinary  and  extraordinary  waves  in  a  mag¬ 
netized  plasma.  Satya  and  Schmidt13  applied  it  to  the  prob¬ 
lems  of  laser  and  Alfven-wave  scattering  from  fluctuations 
in  a  plasma. 

In  the  next  section  we  shall  describe  the  formalism  used 
to  derive  the  general  dispersion  relation  in  the  presence  of 
fluctuations.  Specialization  to  the  case  of  a  whistler  will  be 
the  topic  of  Sec.  III.  Here  we  shall  define  an  attenuation 
length  for  a  whistler  wave,  which  can  be  interpreted  as  a 
distance  over  which  the  coherent  form  of  the  wave  energy  is 
lost.  This  wave  energy  then  reappears  as  an  incoherent  whis¬ 
tler  energy.  Thus  another  way  of  describing  what  is  occur¬ 
ring  is  to  say  that  the  coherence  of  the  beam  is  gradually 
being  eaten  away  by  the  fluctuations.  Next,  in  Sec.  IV,  we 
shall  evaluate  the  attenuation  length  of  the  coherence  for  the 
typical  parameters  of  tokamak  and  space  plasmas.  Finally, 
in  Sec.  V,  we  summarize  the  main  conclusions  and  also  make 
some  pertinent  remarks. 

II.  FORMULATION  AND  BASIC  EQUATIONS 

First  we  shall  state  the  basic  assumptions  involved:  (i) 
Consider  a  cold,  collisionless  and  homogeneous  plasma  im¬ 
mersed  in  a  steady  magnetic  field  directed  along  the  z  axis, 
(ii)  In  order  that  the  perturbation  expansion  be  applicable, 
we  shall  assume  that  the  fluctuations  are  small  over  a  uni¬ 
form  background  plasma,  (iii)  We  shall  treat  the  fluctuations 
as  quasistatic  so  that  the  frequency  characterizing  the  fluctu¬ 
ations  may  be  taken  to  be  much  smaller  than  that  of  the 
whistler  wave,  (iv)  The  fluctuation  spectrum  is  chosen  to  be 
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of  a  Gaussian  shape  for  analytic  simplicity.  This  choice  of  a 
Gaussian  shape  is  further  motivated  by  the  fact  that  the  ob¬ 
served  spectrum  can  be  closely  approximated  by  a  Gaus¬ 
sian.13  (v)  In  keeping  with  the  observations,2,3  we  shall  also 
assume  that  the  fluctuation  spectrum  is  isotropic  in  a  plane 
perpendicular  to  the  direction  of  the  external  magnetic  field 
(i.e.,  the  z  axis).  Further,  since  the  fluctuations  are  mostly 
field-aligned,  i.e.,  the  density  distribution  is  stratified  along 
the  ambient  magnetic  field,  we  shall  ignore  their  effect  on  the 
1  component  of  the  wave  refractive  index  («, ).  We  shall  also 
regard  the  spectrum  as  a  function  of x  and  >>  only  and  take  the 
correlation  length  LT  to  be  the  same  in  these  directions,  (vi) 
We  shall  consider  the  whistler  frequency  <a0  to  be  in  the  low¬ 
er-hybrid  range,  i-e.,  (1,  with  tla  denoting  the  cy¬ 

clotron  frequency  of  species  a.  This  is  the  frequency  range  of 
interest  in  the  rf  heating  of  thermonuclear  plasmas.9  Whis¬ 
tlers  in  this  frequency  range  also  occur  in  space  plasmas. 14 

Assuming  time  harmonic  solutions,  the  Maxwell  equa¬ 
tions  can  be  written  in  the  form, 

VXVXE  =  K-E,  (1) 


where  the  spatial  variables  have  been  normalized  in  units  of 
(Dq/c  and  K  is  the  cold  plasma  dielectric  tensor  given  by 


K  = 


X 

X 

o 


(2) 


where  the  elements  X ,  X ,  and  K  „ 

are  defined  as  follows: 

X  ®  1  —  +  <o\,/Q  i. 

(3a) 

X=»J./®(A. 

(3b) 

X  ss 

(3c) 

Here  to^  is  the  plasma  frequency  of  species  a. 

If  there  were  no  fluctuations,  then  we  know  that  Eq.  (1) 
would  just  provide  us  with  the  usual  linear  description  of  the 
waves.  Because  of  the  presence  of  the  density  fluctuations, 
however,  the  density-dependent  parts  of  the  cold  plasma  di¬ 
electric  tensor  K  are  modified.  Then  one  can  expresses  Eq.  ( 1 ) 
in  the  form  of  a  general  stochastic  equation  {i.e.,  a  partial 
differential  equation  with  random  coefficients), 

0.  (4) 


where  is  the  stochastic  differential  operator  (matrix),  and 
£  is  a  random  wave  field.  It  is  convenient  to  separate  out  the 
part  due  to  fluctuations  from  the  general  Maxwell  operator 
if  in  the  basic  equation  (4).  This  is  done  simply  by  writing 
the  dielectric  tensor  demenu  formally  as, 


Xo  +  (Xo  -  lHSn/AW. 

(5) 

Xoll+fin/Ar,,), 

(6) 

*go(l  +  Sn/ Nq), 

(7) 

where  the  quantity  {6n/N0)  arises  due  to  the  density  fluctu¬ 
ations  and  the  quantitites  with  zero  subscript  are  the  unper¬ 
turbed  parts  evaluated  at  the  ambient  plasma  density. 

With  these  developments,  we  can  split  Eq.  (4)  as  follows: 

i'£-(I  +  K)£«l{l+£.-,K)£.0.  (8) 

Here,  L  represents  the  nonrandom  part  of  the  operator, 
while  V  stands  for  the  fluctuating  part.  The  above  equation 


tacitly  assumes  that  the  inverse  of  L  exists.  The  justification 
of  this  assumption  comes  from  the  consideration  that  the 
fluctuations  alter  the  wave  refractive  index,  so  that  in  their 
presence  the  determinant  of  L  is  nonzero.  Further,  the  fluc¬ 
tuating  part  V  can  be  expressed  as 

V  —  {Sn/N0]M  (9) 

and  M  is  the  constant  matrix, 

/-X 0-1)  X o  o  \ 

M-i  -Xo  -(Xo-1)  0  1.  (10) 

\  0  0  -Kj 

Basically  Keller’s  method  involves  the  derivation  of  an 
equation  for  the  averaged  field  in  the  presence  of  fluctu¬ 
ations."  To  derive  such  an  equation,  begin  with  Eq.  (4)  and 
write  the  stochastic  operator  Jif  as  follows: 

=  [Jf-')"'  =  |(1  +L  ~'Vr'L  (11) 

Next,  consider  the  binomial  expansion  for  ( 1  +  L  ~ 1  V)~ ', 


(1  +L-'V)~'=  £(- L-'V)\  (12) 

H  •»  0 

This  requires  that  \L  ~ 1 V  |  <  1.  Thus,  one  must  estimate  the 
magnitude  of  this  last  quantity  to  ensure  that  the  binomial 
expansion  as  such  will  be  applicable  in  a  given  problem.  For 
now,  it  will  be  assumed  that  \L  ~ 1 V  |  is  sufficiently  small  to 
permit  dropping  of  terms  up  to  and  higher  than  ( L  ~'Vf. 
This  will  be  a  posteriori  justified  later  in  this  paper. 

Next,  rewrite  Eq.  (4)  using  Eqs.  (11)  and  (12)  and  consid¬ 
er  the  ensemble  average  of  the  resulting  equation,  i.e., 

RS'- 

- IL  +  (V)-(VL-'V)  +  (V)L-'(V > 

+  0(L~'V?](E)=  0,  (13) 

where  the  angle  brackets  indicate  an  ensemble  average.  If 
{V)  =  0,  and  the  terms  of  order  (L  ~ 1  Vf  and  higher  are  ig¬ 
nored  in  the  above  equation,  then  one  obtains  Keller's  equa¬ 
tion, 

L  (E)  *=  (VL~'V)(E).  (14) 

The  inverse  operator  £  “ ',  or  the  Green’s  function,  is  defined 
such  that, 


£<7(r,r’)«  W(r-  y’), 


(15) 


where  G(r,r')  is  the  usual  Green’s  function  or  in  this  case  a 
Green’s  matrix.  We  may  also  write  down  the  following 
Fourier  representation  for  Green’s  matrix  G(r,r'),  obtained 
by  simply  taking  the  Fourier  transform  of  Eq.  (IS): 


O(M0 


l£(n  4-  a')£  (n) 

I*  Ml  ’ 


(16) 


where,  £  (n)  represents  the  transposed  matrix  of  cofactors  of 
L  and  \L  fn)|  is  the  determinant  of  L.  In  terms  of  G,  L  ~ '  can 
be  expressed  as  an  integral  operator. 


JdVG(r,r')K(r'). 


(17) 


Consequently,  Eq.  (14)  becomes 


L  (E)  -  (JrfV  V (r)G (r,r'| V{r')(E (r')> (18) 
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From  Eq.  (18)  one  can  readily  derive  the  dispersion  rela¬ 
tion  using  the  standard  Fourier  technique."  The  result  can 
be  expressed  as  follows: 


detfl,., - 

\  (2ir),/2J  |£(n')|  7 


(19) 


In  Eq.  (19),  we  have  defined  a  two-point  correlation  function 
a\n  —  «')  such  that 


_<gn(r)  jn(rj>  = 

AT0  N0 

i.e.,  the  medium  is  statistically  homogeneous.  We  could 
solve  Eq.  (19)  for  different  forms  of  the  correlation  spectrum 
such  as  Gaussian,  exponential,  Loren tzi an,  etc.  As  men¬ 
tioned  in  assumption  (iv)  at  the  beginning  of  this  section,  we 
shall  illustrate  the  calculations  for  a  Gaussian  spectrum. 
Thus  take 


o(n  -  n')  =  r*72  exp[|nt  -  n[  1 *i2/2]6{nz  -  n’z),  (20) 

where  =  kLc/a> 0  is  the  perpendicular  component  of  the 
propagation  vector  and 7  =»  <uo/LT/c  is  the  normalized  cor¬ 
relation  length.  Further,  e  =  ((<5rt/n0)2) 1/2  is  the  measure  of 
the  strength  of  the  fluctuations.  Note  that  if  e  =  0,  i.e.,  if 
there  were  no  fluctuations,  the  solution  to  Eq.  (19)  would 
yield  the  usual  dispersion  relations  of  the  characteristic 
modes. 


III.  ANALYSIS  OF  THE  GENERAL  DISPERSION 
EQUATION 

First,  w  remark  that  Eq.  ( 19)  contains  two  cold  plasma 
modes.  One  of  these  is  the  “slow”  or  lower-hybrid  mode  and 
the  other  one  is  the  “fast”  or  the  whistler  mode.  Our  object  is 
to  focus  on  the  latter  mode.  To  this  end,  we  need  only  to 
recall  that  the  whistler  mode  polarization  is  such  that  |£z  j 
<\E,  |  or  (£,  |,  i.e.,  we  can  neglect  |£z  |  compared  with  either 

_ I 


(£,  |  or  |£j,  |  in  the  frequency  range  of  interest.11  Then  we  can 
approximate  the  3  X  3  matrix  L  to  the  following  2x2  matrix 

Klq  —n\—n\  sin2  d  —  iKM  +  rtf  sin  6  cos  4> 
iKM  +  rt\  sin  d  cos  d  Kl0  —  n2  —  n2  cos2  d 

121) 

The  corresponding  matrix  M  is 


/  -  (Kl0  -  1)  K*  \ 
\-iKM  -\Kl0  -  1)/ 


(22) 


We  note  from  Eq.  (21)  that  the  whistler  dispersion  relation 
follows  from  the  solution  to  \L  (n)|  =  0,  i.e., 

«20^£io/(«2-£i0)-(n2-Ari0) 

■sK\o/n\ -n\,  if  n\>Kl0.  (23) 

To  evaluate  the  integral  in  Eq.  ( 1 9),  we  choose  a  coordi¬ 
nate  system  such  that  n,  =  nLx,  and  let  d  be  the  scattering 
angle.  Then  the  exponent  in  Eq.  (20)  can  be  written, 

|»t  —  nj  |  =  n\  +  n[2  —  lnLn[  cosd.  (24) 

Also,  dri  —  drt'x  dn'ydn'z  =  n[  dn[  drt '  d<j>.  With  these  preli¬ 
minaries  we  can  evaluate  the  n[  -integral  in  Eq.  ( 19)  by  noting 
that  the  integrand  is  singular  at  the  roots  of  | L  |n')|  =  0,  on 
the  real  axis.  For  a  forward  propagating  incident  wave,  the 
appropriate  contour  passes  an  infinitesimal  distance  below 
the  pole.  Then  the  imaginary  contribution  can  be  simply 
found  by  using  the  Cauchy  principal  value  formula.  The  in¬ 
tegral  over  n’  is  trivial  due  to  the  delta  function.  Finally,  for 
the  angle  integration  we  use  the  integral  representation  of 
the  modified  Bessel  functions.15  After  evaluating  the  inte¬ 
grals  for  each  matrix  element  in  Eq.  (19)  and  some  algebra, 
we  can  write  the  expressions  for  the  real  and  imaginary  parts 
of  the  complex  refractive  index  a  L  in  the  following  compact 
forms: 


ReK,=r° - 


-Im(ni)s^/,2  + 


[2£,p/|a(n2  —  Kl0)  +  /|)(£ —  /  i;  +  f 1 7;;) ] 

(n2 


)(,+ra) 


(25) 

(26) 


In  Eqs.  (25)  and  (26)  the  quantities  t,,  arise  due  to  the  fluctu¬ 
ations  and  are  defined  in  the  following  relations: 

(27) 

/,,=  -/I( 

=  A  -  l)[*io  +  (a*  -  ATgjj K*  -  1)] 

2  (nj  —  fCl0  )2 

(28) 

I„-A{a'FM-m*\l  (29) 

where 


r 


£0(s)  =  I0{r)e  *\ 

£,(*)  = /,(s2)s_2*_'2, 

A  =  (,r/2)"  V72), 

and  /0  and  /,  are  the  modified  Bessel  functions  order  zero 
and  one,  respectively.  Their  argument  is  s2  =  n2l0l 2,  nl0  be¬ 
ing  the  unperturbed  perpendicular  refractive  index  given  by 
the  linear  dispersion  relation  (23),  and  7  was  defined  just  be¬ 
low  Eq.  (20).  We  have  also  defined 


-*>j-l) 

(n|-£l0)2 
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a  = - ; - ; - ,  (Jl) 

(«r  -Kio) 

' - 5T^F - '•' 

Naturally  one  may  note  from  the  relations  (25)  and  (26)  that 
in  the  absence  of  fluctuations,  the  7,/s  vanish  giving 
Re(n1 )  =  n10  and  Im(nt )  =  0  as  expected. 

Equation  (23)  expresses  the  change  in  the  real  part  of  the 
refractive  index  due  to  the  fluctuations.  This  in  turn  means 
that  the  wave  velocity  experiences  only  second-order 
changes  due  to  density  perturbations.  On  the  other  hand,  Eq. 
(26)  says  that  the  primary  wave  is  being  attenuated  by  the 
presence  of  the  fluctuations  as  the  wave  energy  scatters  off 
them  and  is  a  first-order  change.  One  can  define  the  attenu¬ 
ation  length  as  follows: 

/„  =c/[2*>0ImK)].  (33) 

This  is  the  distance  over  which  the  coherent  energy  from  the 
primary  wave  is  lost  to  the  fluctuating  components. 

Before  discussing  applications  of  the  basic  results  of  this 
section,  namely  Eqs.  (23)  and  (26),  let  us  consider  two  limit¬ 
ing  forms  of  these  equations.  These  correspond  to  the  small 
and  large  argument  limits  of  the  modified  Bessel  functions. 13 

(i)  i2<7,  or  the  long  wavelength  limit-.  In  this  case  we  can 
use  the  small  argument  expansions  of  the  Bessel  functions. ,s 
We  thus  have,  F0(s)  ail  and /■,($)  ai  1/2.  Weshallfind  that  for 
most  applications  we  have  in  mind,  the  dielectric  tensor  ele¬ 
ment  J^x0>l  and  that  >/CL0 .  Further,  it  can  be  shown 
that  n]  *  Km  cos  9^  where  ff0  is  the  propagation  angle  rela¬ 
tive  to  the  background  magnetic  field.  Using  this  fact  and 
retaining  only  the  dominant  contributions,  these  arguments 
simplify  Eqs.  (23)  and  (26)  as  follows: 

Re(nJ )  =  it  Jo  +  O  (4  2),  (34) 

Im(ni)  =  11  +.±^2  **  +  *„)  +  0{A  2).  (35) 

2  cos4  90 

(//)  s2>  l,  or  the  geometric  optics  ( short  wavelength  )  limit-. 
In  this  case,  the  random  inhomogeneity  scale  length  is  large 
compared  to  the  whistler  wavelength.  Thus  we  can  use  the 
asymptotic  forms  of  the  Bessel  functions, 13 

where  v  =  0,1.  Then  we  have 

FM~{2, ra2)1'2’  ^’~7(2^P 

In  this  case  we  find 


TABLE  I.  Typical  plasma  parameters. 


Parameters 

Tokamak 

Ionosphere 
if  layerl 

PLasmasphere 
(4  earth  radii  | 

Whistler 
frequency 
u„rad  sec  1 

5x  lO4 

5  xIO4 

5x  103 

Electron 

density 

.Voicm”  *1 

2X10° 

10“ 

10* 

Magnetic 

field. 

B,  (Gaussl 

2:<  I04 

4.8  x  10-' 

5X10-' 

Electron 

temp. 

TJtV) 

101 

1.5x10-' 

1 

Ion 

temp. 

TA  eV} 

200 

3  x  10  : 

1 

s2  is  such  that  it  satisfies  the  condition  L  ~ 1 V :  <  1 .  More  on 
this  will  be  said  later. 

We  shall  now  discuss  the  implications  of  these  results  in 
the  next  section. 

IV.  APPLICATIONS 

To  discuss  the  implications  of  the  results  of  the  last  sec¬ 
tion,  let  us  begin  by  listing  the  typical  plasma  parameters  of 
interest  in  Table  I.  A  glance  at  Eqs.  (25)  and  (26)  reveal  that 
the  knowledge  of  the  dielectric  tensor  elements  KM  and  K.0 
is  also  required.  The  estimates  of  these  in  the  frequency 
range  of  interest  are  given  in  Table  II. 

For  the  subsequent  discussion  it  will  be  helpful  to  recall 
from  the  whistler  propagation  theory8  that  in  the  chosen 
frequency  regime,  the  whistler  propagation  angles  can  have 
a  wide  range  of  values  with  the  upper  bound  determined  by 
the  condition  that  the  refractive  index  remain  real.  On  the 
other  hand,  the  group  velocity  (ray)  direction  remains  within 
a  maximum  of  19.5  deg  relative  to  the  background  magnetic 
field.  In  what  follows  we  shall  examine  the  behavior  of  the 
attenuation  length  as  a  function  of  the  propagation  angle  of  a 
whistler  wave. 

The  next  task  is  to  have  the  estimates  of  the  fluctuation 
level  e  and  the  correlation  length  LT.  For  this  purpose  we 
need  to  invoke  a  suitable  model.  We  shall  choose  a  model 
deduced  from  observations2  '  according  to  which  the  fluctu¬ 
ations  are  caused  by  drift-wave  type  turbulence.  Drift  waves 
are  the  characteristic  modes  of  an  inhomogeneous  plasma. 


InsilA  /[hrs1)11^  +  0  /s1)],  (36) 

/  i  ~  [Kgo  cos  0O  —  1) 

7,2  21  “  ^W*2)1'2  COS2  00  ’ 

/n-fii/M'^Ua'-tf/i2).  (38) 

With  these,  Eqs.  (23)  and  (26)  become 

Re(«i)-«io  +0[A\  (39) 

Im(n2)  *  [/Itf  lo/Utrr2)1'2  cos2  0O)  +  0(4  2).  (40) 

It  should  be  kept  in  mind  that  this  limit  is  meaningful  only  if 
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The  extremely  low-frequency  (ELF)  electric  field  fluctu¬ 
ations  observed  near  the  plasmapause  (boundary  of  the  plas- 


TABLE  It.  Estimates  of  dielectric  tensor  elements  in  LH  range 


Ionosphere 

Plasmasphere 

Element 

Tokamak 

If  layerl 

(4  earth  radii  I 

36 

7.5x10' 

7.2  x  10' 

*t0 

0.12 

6.3 

35 
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TABLE  HI.  Estimates  of  t  sad  Lr  for  drift-wave  model. 


Parameter 

Tokamak 

Ionosphere 
[F  layer) 

Plaemaaphere 
(4  earth  radii) 

L,  cm 

10 

10’ 

6.4  x  10' 

Lr  cm 

7.2x10-* 

1.6x10* 

2  x  10* 

e 

7.3x10-’ 

1.6x10-’ 

3.1  XlO'4 

masphere)  region  can  be  explained  in  terms  of  the  drift  insta¬ 
bility.16  In  the  equatorial  F  region  of  the  ionosphere,  drift 
waves  presumably  dominate  for  wavelengths  A  S  100  m.1 

Thus  adopting  this  drift-wave  model  we  can  assume2,3 
the  following  empirical  scaling: 

LT~p„  (41) 

e~L,/Ln,  (42) 

where  pi  —  u„/Fl,  (v„  is  the  ion  thermal  speed)  is  the  ion 
gyroradius  and  L„  is  the  scale  length  of  the  regular  density 
gradient.  For  the  present  discussion,  we  have  L„  >LT.  Table 
III  presents  estimates  of  e  and  L  T  for  the  chosen  parameters. 

A.  Case  1.  Tokamak  plasma 

First  we  note  that  in  order  for  the  plasma  interior  to  be 
accessible  to  the  externally  launched  rf  wave,  the  parallel- 
refractive  index  is  limited  to  a  rather  narrow  spectrum.9 
Typically,  1.3  <  n,  <  3.  Now  for  the  whistler  wave,  we  have 
nz  —  Kxo  cos  ff0  so  that  the  accessible  nx  values  correspond 
to  the  propagation  angle  range  of  75.6-87.3  deg  for  the  KM 
value  appearing  in  Table  II.  Estimating  the  perpendicular 
refractive  index  #il0  from  the  whistler  dispersion  relation  for 
the  allowed  nz  values  and  taking  the  value  for  the  correlation 
length  Lr  from  Table  III,  we  find  thatr2  =  (n107)2<l.  Thus, 
we  use  the  appropriate  limiting  forms  of  the  real  and  imagi¬ 
nary  parts  of  the  refractive  index  given  by  Eqs.  (34)  and  (35), 
respectively.  To  the  leading  order,  we  see  that  the  real  part  of 
the  refractive  index  is  unchanged,  while  the  imaginary  part 
can  be  used  in  the  definition  of  the  attenuation  length  [Eq. 
(33)]  to  get 

yl  0°  )  (43) 

Equation  (43)  makes  the  0o  dependence  of  the  attenuation 
length  explicit.  Figure  1  shows  the  plot  of  the  normalized 
attenuation  length  (4,/a),  a  being  the  minor  radius  of  the 
tokamak.  We  notice  from  this  plot  that  the  ratio  (/„  /a)  drops 
as  a  increases  in  the  chosen  range.  This  indicates  increased 
scattering  from  the  density  fluctuations  and  the  consequent 
loss  of  coherence  of  the  primary  whistler  wave. 

8.  Cbm  2.  Ionospheric  (flayer)  plasma 

Consider  a  typical  whistler  propagating  in  the  F  region 
of  the  ionosphere.  Calculating  nL0  from  the  linear  dispersion 
relation  for  different  values  of  the  propagation  angles  80  and 
picking  the  Lr  value  from  Table  III  under  the  ionosphere 
column,  we  again  find  that  j2  <  1 ,  thus  allowing  the  use  of  the 
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FIG.  I.  Normalized  attenuation  length  versus  angle  I  tokamak  plasmal. 

approximate  Eq.  (43).  In  Fig.  2  we  have  plotted  the  normal¬ 
ized  attenuation  length  (7a  =  as  a  function  of  the 

propagation  angle  60. 

We  notice  from  Fig.  2  that  /„  decreases  with  the  in¬ 
crease  of  60.  For  the  model  being  considered,  /„  is  of  order 
103  km. 

C.  Cass  3.  The  plasmapaus#  region 

The  drift  waves  presumably  dominate  in  this  region  lo¬ 
cated  at  about  4  earth  radii.  Thus  we  can  again  apply  the 
fluctuation  model  being  considered.  Using  the  LT  values 
from  Table  III  and  knowing  nl0  from  the  dispersion  rela¬ 
tion,  we  find  that  the  long  wavelength  limit  is  again  ade¬ 
quate.  Thus  making  use  of  Eq.  (43),  we  find  that  the  behavior 
of  7a  vs  0O  is  qualitatively  the  same  as  in  ionosphere,  but  here 
/„  is  of  order  10*  km.  Comparison  of  this  value  with  the 
correlation  length  and  the  scale  length  of  the  regular  density 
gradient  given  in  Table  III  shows  that  is  larger  than  both 
these  lengths.  Within  the  fluctuation  model  being  consid¬ 
ered,  we  would  expect  then  that  the  loss  of  coherence  of  a 
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FIG.  2.  Attenuation  length  venue  angle  i ionosphere) 
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whistler  wave  due  to  scattering  from  fluctuations  in  the  plas- 
mapause  region  to  be  insignificant. 

V.  CONCLUSIONS  AND  REMARKS 

In  this  paper  we  have  considered  the  scattering  of  a 
whistler  wave  from  random  density  fluctuations  in  a  magne¬ 
tized  plasma.  We  derived  a  modified  dispersion  relation  of 
the  wave  by  using  Keller’s  method  assuming  the  randomness 
to  be  small  compared  to  the  nonrandomness.  Although  we 
discussed  only  the  case  of  a  Gaussian  correlation  spectrum, 
other  forms  could  be  readily  included. 

The  expression  for  the  attenuation  length  was  simpli¬ 
fied  and  was  examined  for  the  scattering  from  short  wave¬ 
length  fluctuations.  The  results,  in  particular,  showed  that 
scattering  becomes  important  as  the  propagation  angle  in¬ 
creases.  We  should  remark  here  that  the  model  used  to  esti¬ 
mate  the  correlation  length  and  the  level  of  fluctuations  is 
not  expected  to  be  always  valid.  In  particular,  quite  large 
fluctuation  levels  (about  30%  or  more)  have  been  reported  in 
the  Alcator  tokamak.2  When  the  fluctuations  become  this 
strong,  the  perturbation  method  used  here  may  of  course  be 
questioned.  However,  these  results  may  still  be  used  to  deve- 
lope  a  qualitative  understanding.  Nevertheless,  for  any  de¬ 
tailed  or  rigorous  analysis,  there  is  definitely  a  need  for  some 
techniques  designed  to  handle  strong  fluctuations. 

Also,  at  the  present  time  it  is  not  clear  how  to  predict,  in 
a  reliable  maimer,  the  fluctuation  levels  in  the  future  reactor 
size  devices.  Consider,  for  example,  a  reactor  size  D-T  plas¬ 
ma  device  with  fluctuation  level  of  1  %,  minor  radius  =  270 
cm,  iV0  =  10M  cm'3,  Ba  =  32  kG,  nj  =  10,  T,  =  8  keV, 
eo0  =  J2c0lh ,  where ailh  =  «jp, ( 1  +  2)~ 1,2  is  the  lower- 

hybrid  frequency.  For  this  case,  r  =  5.6  and  reading  the 
Bessel  function  values  from  the  table, 15  we  find  that  the  ratio 
(/„  /a)  is  of  order  10"  \  indicating  severe  attenuation  due  to 
scattering  before  the  wave  has  a  chance  to  penetrate  even  one 
percent  of  the  distance  to  the  center.  This  essentially  stems 
from  the  large  size  of  the  device.  We  also  remark  that  fluctu¬ 
ations  of  other  kinds,  such  as  those  in  the  magnetic  field,  also 
occur  in  tokamalcs  and  these  need  to  be  included  in  a  more 
complete  theory.  But  these  calculations  do  demonstrate  the 
importance  of  scattering  from  random  density  fluctuations. 

As  one  can  note,  the  attenuation  length  is  obtained  from 
the  modified  dispersion  relation  of  the  wave  in  the  presence 
of  the  fluctuations  and  it  does  not  carry  any  information  on 
the  evolution  of  the  wave  field  itself.  To  gain  an  insight  into 
this  question  of  how  the  field  actually  evolves  in  the  presence 
of  fluctuations,  one  must  consider  the  evolution  equation  for 
the  field.  Start  for  this  purpose  from  the  equation  for  the 
averaged  field  given  in  Eq.  (14)  and  perform  the  standard 
multiple  scale  analysis  on  this  equation.17  As  a  result  one 
obtains  the  following  linear  equation: 

id,*+P%4  +  Qt  =  0.  (44) 

This  equation  assumes  that  the  parallel  refractive  index  of 
the  wave  is  not  changed  due  to  the  fluctuations.  Further,  the 
coefficient  P  stands  for  the  dispersion,  and  is  defined  as 


while  Q  is  the  pure  imaginary  loss  term  due  to  attenuation. 
Equation  (44)  describes  the  linear  spatial  evolution  of  the 
wave  field  taking  account  of  the  attenuation  due  to  scatter¬ 
ing.  This  equation  is  reminiscent  of  the  neutron  diffusion 
equation  in  a  nuclear  reactor18  where  it  describes  the  slowing 
down  and  diffusion  of  neutrons  in  a  moderator  with  the  loss 
term  signifying  the  absorption.  In  the  present  discussion  no 
diffusion  occurs  in  the  limit  considered,  but  only  pure  at¬ 
tenuation  of  the  primary  whistler  occurs  as  it  scatters  off  the 
density  fluctuations. 

From  what  has  been  said  thus  far.  we  only  know  that 
some  of  the  coherent  wave  energy  lost  due  to  scattering  gets 
transformed  into  the  incoherent  form.  It  might  be  of  interest 
to  ask  how  does  one  calculate  this  incoherent  part  of  the 
intensity.  For  this  purpose  one  would  have  to  reconsider  the 
basic  equations  (4)  and  ( 14),  noting  that  E  in  14)  is  the  total 
field  and  (E )  occurringinEq.il  4)  is  the  same  as  the  coherent 
field.  Then  writing  E  =  E,  +  Ec,  and  taking  the  difference 
of  these  equations,  one  gets  the  equation  for  the  incoherent 
field.  It  is  this  latter  equation  that  one  must  solve  to  find  the 
incoherent  part  of  the  field. ' 1  Note  that  in  the  model  being 
considered  the  energy  in  the  fluctuating  components  would 
remain  in  the  wave  field  only.  Its  transfer  to  the  particles 
could  be  predicted  not  from  the  present  model  but  from  a 
model  that  would  include  the  appropriate  wave-particle 
coupling.  Thus  within  the  present  model,  imagining  the  pri¬ 
mary  whistler  as  a  radiating  antenna,  one  could  say  that  as 
the  wave  propagates  in  the  presence  of  fluctuations,  it  would 
radiate  incoherent  waves  at  the  same  frequency  but  different 
harmonics  of  the  wavenumbers.  This  radiated  wave  energy 
would  be  removed  from  the  main  whistler  and  would  per¬ 
vade  the  rest  of  the  plasma  in  a  turbulent  form. 

Finally,  consider  the  question  of  the  range  of  validity  of 
Keller’s  method.  As  has  been  previously  noted,  one  must 
estimate  the  magnitude  of  the  quantity  j  L  ~ '  V  \  for  this  pur¬ 
pose.  Now  L  ~ 1  is  some  function  of  the  wave  refractive  index 
modified  by  the  fluctuations  (i.e.,  if  n  is  the  effective  index  of 
refraction  with  n  =  n°  +  n\  where  n°  is  the  unperturbed  val¬ 
ue  and  n’  is  the  small  perturbed  part  introduced  by  the  fluc¬ 
tuations).  Then  expanding  L  ~  \  one  finds  to  the  leading  or¬ 
der  that  L  ~ 1  is  of  order  (Im  nL/nL ).  Now  the  fluctuating  part 
V  is  of  order  e,  where  e  denotes  the  strength  of  the  fluctu¬ 
ations.  Thus  one  has  \L  ~ 1 V  \  to  be  of  order  (Im  nL/nL]e.  For 
the  whistler  case,  from  the  matrix  (21)  one  can  show  that 
|  L~'V\  goes  like  (se^/f  ^ ,  so  that  the  condition  L  ~ 1 V  |  <  1 
becomes 

sKUf1^.  (45) 

This  sets  a  limit  to  the  value  of  the  parameter  s  up  to  which 
the  binomial  expansion  (12)  could  be  safely  applied.  Taking 
the  typical  KM  values  listed  in  Table  II  and  the  values  of  the 
fluctuation  level  e  as  given  in  Table  HI,  we  find  that  Keller’s 
method  is  applicable  for  s  values  ranging  below  about  10- 
100. 
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